Interplay between interference and Coulomb interaction in the ferromagnetic 
Anderson model with applied magnetic field 
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We study the competition between interference due to multiple single-particle paths and Coulomb 
interaction in a simple model of an Anderson-like impurity with local-magnetic-field-induced level 
splitting coupled to ferromagnetic leads. The model along with its potential experimental relevance 
in the field of spintronics serves as a nontrivial benchmark system where various quantum transport 
approaches can be tested and compared. We present results for the linear conductance obtained 
by a spin-dependent implementation of the density matrix renormalization group scheme which are 
compared with a mean-field solution as well as a seemingly more advanced Hubbard-I approximation. 
We explain why mean-field yields nearly perfect results, while the more sophisticated Ifubbard-I 
approach fails, even at a purely conceptual level since it breaks hermiticity of the related density 
matrix. Furthermore, we study finite bias transport through the impurity by the mean-field approach 
and recently developed higher-order density matrix equations. We find that the mean-field solution 
fails to describe the plausible results of the higher-order density matrix approach both quantitatively 
and qualitatively as it does not capture some essential features of the current-voltage characteristics 
such as negative differential conductance. 

PACS numbers: 72.25.-b,85.75.-d,73.23.Hk,73.63.-b 



INTRODUCTION 



When electrons pass through a mesoscopic region, the 
superposition of several different single-particle transport 
paths can lead to interference, as, e.g., in an Aharonov- 
Bohm geometry with quantum dots embedded in the 
armsi^ As the size of the mesoscopic region dimin- 
ishes, many-particle effects such as Coulomb blockade be- 
come increasingly important.'^ This may change the am- 
pUtudes of the competing transport paths and thereby 
alter the interference effect. Eventually, for sufficiently 
strong many-body interaction, the single-particle-paths 
picture breaks down and such systems should be treated 
using a true many-body formalism. This problem of 
the interplay between interference of several competing 
paths and many-body interaction has recently attracted 
a lot of attention theoretically in the general quantum- 
transport contexi^^^^^^^^S^^^^^^^^^^^^^ as well as 
from more specific points of view such as the molec- 
ular electronics^^'21^'23^ spintronics ^^^^s^ 
or even full counting statistics'^" and superconducting 
transport.— 

In the case of spintronics, the interference can be 
achieved without necessity of a multiply-connected or- 
bital geometry due to the possibility of superposition 
of different purely spin amplitudes with the help of ci- 
ther non-coUinearly magnetized lead a^^i^^'^^ or an ad- 
ditional spin-level splitting non-coUinear with the lead 




FIG. 1: Sketch of our model, where the magnetic field B in 
the central region is tilted by an angle (j) with respect to the 
magnetization of the contacts. 



magnetization ] ^^1'^^ Since experiments with strongly- 
interacting quantum dots and ferromagnetic contacts 
have recently been successfully performed^i^i^^ the spin 
interference effects proposed in previous work^ and fur- 
ther elaborated here may be within experimental reach. 

Previously, some of us considered a model^ consisting 
of a spin- 1/2 level coupled to ferromagnetic leads with 
the magnetizations being either parallel or antiparallel. 
In addition, a magnetic field non-coUinear with the spin 
direction of the leads was applied (see Fig.[T|).®'^ For this 
Ferromagnetic Anderson model with an applied magnetic 
field B (from now on nicknamed the FAB model) the 
linear conductance was obtained in two different regimes: 
without interactions on the dot, and in the cotunneling 
regime. For the non-interacting case with fully polarized 




FIG. 2: The linear conductance of Eq. ((2|, which shows the 
anti-resonances for non-interacting electrons (full line, (7 = 0) 
and the spin valve behavior for strongly interacting electrons 
(dashed line, U = lOOB) as a function of the angle between 
the magnetizations of the leads and the applied magnetic field. 
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FIG. 3: Schematic energy spectrum in the linear conduc- 
tance regime for the non-interacting case. The bare resonant 
level is split due to magnetic field, and the angle between the 
magnetizations of the leads and the applied magnetic field is 
denoted by (j). The widths of the two resonances depend on (j) 
as cos^((^/2) and sin^(0/2), respectively. 



leads, zero temperature, the bare level on resonance and a 
parallel lead configuration, the linear conductance can be 
calculated, e.g., using non-equiUbrium Green functions 
(NEGF), and the exact result ia^ {h is Planck's constant) 



R 



h l + cos2 0[(ro -f r?j)/2B]' 



(1) 



with being the coupling to the leads with a = L,R, 
B the magnetic field times the magnetic moment of the 
level, and the angle between the magnetization in the 
leads and the applied magnetic field. The conductance 
shows anti-resonance at angles (f> = 7t/2 and 4> = 3tt/2 
due to destructive interference, see Fig. 

Under the same conditions as stated above, the linear 
conductance can be obtained in the cotunneling regime 
([r^ -|- T^]/B <^ 1) even in the presence of an on-site 
Coulomb interaction U, see Eq. ([5]), by applying a scat- 
tering formalism:'^^ 



6"= 



cos2(0/2) sin2(0/2) 



-B 



B + U 



(2) 



In this regime, the conductance shows a cross-over from 
the behavior with anti- resonances at (j) = 7r/2 for the 
non-interacting case, to a spin-valve effect for U —^ oo 
with C^ot"" cx cos''((/)/2). That is the anti-resonances 
around (p — Tr/2 disappear and the conductance vanishes 
for (j) = TT instead (see Fig. [2]). 

A simple physical picture for the situation described 
above is: In a basis where the Hamiltonian for the iso- 
lated dot is diagonal, the bare dot level energy is split by 
the magnetic field, and for non-interacting electrons the 
density of states has peaks at the two single-particle en- 
ergies at ^B, see Fig.[3l The widths of the two peaks de- 
pend on the angle ((>, and for fully polarized leads they are 
proportional to cos'^{<j)/2) or sin^((/)/2), respectively. For 
(/) = and (j) = TT one of the peaks is infinitely narrow and 



electrons can only pass through the other level, whereas 
for = 7r/2 the peaks are equally wide resulting in the 
sharp anti-resonances due to interference. So the angu- 
lar dependence of the conductance can be understood as 
interference through non-degenerate levels, which have 
widths depending on the angle between the magnetiza- 
tions of the leads and the applied magnetic field. For 
a large on-site Coulomb interaction some weight of the 
density of states is moved away from the single-particle 
energies and away from the Fermi level, thereby destroy- 
ing the anti-resonances. 

The qualitative difference between the interacting and 
non-interacting regime is important, as it shows a very 
crucial feature in the transport through mesoscopic sys- 
tems, namely that it is generally not the single-electron 
transport paths which determine the transport, but 
rather many-electron processes. Besides from leading to 
interesting physical effects, it also puts strong demands 
on the theoretical transport formalism applied to such 
systems, as it should be able to handle both the coher- 
ence and the interactions. It also has to be applicable for 
sufhciently low temperature, because otherwise thermal 
fluctuations will wash out the interference effect. That 
makes our model an excellent benchmark for transport 
formalisms. Including interactions, if they are sufficiently 
strong, is a challenge in the standard NEGF formalism 
where all single-particle effects including the interference 
are captured exactly. On the other hand, the density ma- 
trix language (generalized master equation; GME) start- 
ing from exact many-body states of the system, thus in- 
cluding the interaction exactly, faces problems when the 
broadening due to the leads comparable with level split- 
ting (leading to interference effects) is to be incorporated. 
Thus, this kind of models poses significant challenges to 
standard transport approaches even outside notoriously 
difhcult strongly-correlated regimes such as the Kondo 
regime. 

Therefore, we use this model for a detailed compar- 
ison study of the performance of different transport 
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FIG. 4: Sketch of the different parameter regimes. For van- 
ishing interaction U, nonequilibrium Green functions pro- 
vide the complete solution of the transport problem, see 
Eq. IT}. For large level splitting, the scattering formalism al- 
lows for a quantitative description of the cotunneling events, 
see Eqs. ([2)). (|AHl . In this work we provide results for the more 
intricate region of moderate level splitting and a finite on-site 
Goulomb interaction. 



formalisms in the potentially problematic and so far 
not addressed regime of broadening comparable to the 
level splitting i? w and arbitrarily strong interac- 
tion U. In particular, we test the results of higher- 
order, i.e. beyond mean-field, decoupling schemes based 
on NEGFiii^e^ and/or many-body-states-based NEGF 
(Hubbard operator NEGF) approaches .■^^'^'^'^'^i^ 
We find that the Hubbard-I approximation in the frame- 
work of NEGF,*°'^^ frequently applied to the Ander- 
son model with or without ferromagnetic leads fiii^^i^i^ 
gives unphysical and even mathematically wrong results 
for the model considered in this paper. This finding raises 
serious questions about the very foundation of the many- 
body-states-based NEGF approaches i^Si^i^ 

As stated above, so far only the non-interacting and 
the cotunneling regime have been considered (see Fig. [4]), 
and for the latter only in linear response. In this paper, 
we calculate the linear conductance at zero temperature 
for arbitrary values of the tunneling coupling, applied 
magnetic field and on-site Coulomb interaction using a 
density matrix renormalization group (DMRG) scheme, 
see Sec. IIIIl Surprisingly, in the linear response regime 
the results obtained using the DMRG scheme can, in 
certain situations, be reproduced using Green functions 
with a mean-field approximation, which is discussed in 
Sec. lIVl The unexpected failure of the Hubbard-I approx- 
imation in the framework of NEGF^-^ is analyzed in 
Sec.|Vl In section |Vl] we extend the calculations beyond 
linear response by applying a generalized master equa- 
tion formalism^^, which works in a basis of many-particle 
states and takes into account higher-order tunneling pro- 
cesses. In Sec. I VIII the failure of the mean-field Green 
function method for finite bias is demonstrated. Finally, 
we conclude on our findings in Sec. IVIIII Appendix |^ 
contains the cotunneling expression for less than full po- 
larization of the leads and off-resonant transport, and 
App. [B] presents details of the mean-field Green function 
calculation. 



II. THE MODEL SYSTEM 

The model Hamiltonian of the quantum dot coupled 
to magnetic leads is 



(3) 



where 



a—L,R,ka 

Here CT =t / i is the spin of the electrons, a denotes the 
left or right electrodes, which are assumed to be polar- 
ized along the z-axis (the spin quantization axis), either 
parallel or anti-parallel. However, in this paper we only 
consider parallel magnetizations of the leads. The quan- 
tum dot is subjected to a magnetic field B, which is tilted 
by an angle (j) with respect to the z-axis and lies within 
the a;z-plane. Note that we neglect the negative sign of 
the electron charge for simplicity. Thus, the energeti- 
cally preferred spin direction is pointing in the direction 
of B throughout this paper. The dot-Hamiltonian reads 
(n„ = dld„) 

Hd =Y.^odUa + Un^n^ - t^sB ■ T^^,dU,, , (5) 



where is the orbital quantum dot energy, B = l/isBj 
represents the magnetic field splitting, r is a vector con- 
taining the Pauli spin-matrices, and U is the on-site 
Coulomb interaction for double occupancy. In a spin ba- 
sis parallel to B, the dot Hamiltonian is diagonalized as 



Hj,=J2iC,-aB)dld, + Un.nr, 



(6) 



where the do- and d^ operators are related by the unitary 
rotation 

-a _ f cos((/)/2) sin((/)/2) \ , , 

^'T- l^^aada, ±t - I _ gin((/)/2) C0s(^/2) i ' ^' ' 

Finally, the tunneling Hamiltonian is 



a.—L,R k(7 



a—L,R kaa 



(8) 



Here we allow for the tunneling matrix element ta^k„ to 
be spin-dependent, because the states in the leads depend 
on the spin direction. Note that there is no spin-flip asso- 
ciated with the tunneling here, i.e., there is no spin-active 
interface, which would require the use of a non-diagonal 
tunneling matrix ta^kaa' ■ Depending on the parameters 
this would correspond to having an angle between the 
lead magnetizations, which would modify the details but 
not the general behavior that we discuss. 
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We define the energy-dependent coupling constants as 
r^e) ^27tJ2 ~ ^o.m) = E ^".-(e). (9) 

and let Pa denote the polarization of the tunneling from 
lead a defined through Ta^ais) = 5 (1 + <^Pa) ^ai^)- No- 
tice that Pa G [—1,1] such that = ±1 corresponds to 
full spin-t / J, polarization and Pa — corresponds to un- 
polarizcd leads. For parallel (antiparallel) polarization of 
the leads the P^'s have the same (opposite) sign. 

In the basis where the dot part of the Hamiltonian is 
diagonal, the coupling of the two dot states, t and |, to 
the lead a is given by a matrix in the spin index, see also 
Eq. IHSl), 



is) 



Ta{e) 



(1 + aPa COS I 

Pa sin (j) 



for a = a' 
for a ^ a' 



(10) 

In the calculations using the DMRG and the density 
matrix technique, we use a polarization of both leads 
less than 1 for technical reasons. The minority spin only 
introduces a smearing of the results discussed for fully 
polarized leads. 



III. LINEAR RESPONSE: DMRG 

A. Tight-binding Hamiltonian 

In order to apply the DMRG method to the model a 
discretized version of the leads must be formulated. The 
simplest choice is to model the leads as one-dimensional 
semi-infinite tight-binding (TB) chains that are dis- 
cretized appropriately. With this choice and denoting the 
hopping matrix element between the resonant level and 
the leads by ta,a (a = L, R), the tight-binding Hamilto- 
nian reads H^^ = J2a=L.R + + ^d, where 



„TB ^ _ V V — f c -I- c 



(11) 

(12) 



a=L,R a 



and where Ho is given in Eq. ([5]). That is in the DMRG 
implementation we work in the lead spin-basis and do 
not use the diagonalized version of the dot part of the 
Hamiltonian. In i?^^, 2D is the bandwidth of the tight- 
binding chain representation of the leads corresponding 
to the hopping amplitude D /2 between the internal sites 
in the chains. 

In order to link the different theoretical approaches ap- 
plied to solve the model, an expression for the effective 
energy-dependent coupling constants between the single 
site and the tight- binding leads, F™ (e), must be estab- 
lished. For the one-dimensional tight-binding model of 



the leads, these are given by^ 

F™ (e) = -2\taA^ Im [3^^,(1, 1, e)] , (13) 

where (7^g.(l,l,e) is the surface component of the re- 
tarded Green function of the semi-infinite left or right 
chain at energy e. The surface of the tight-binding chain 
is the first site, and the Green function reads^ 



,(l,l,z) = 2- 



Vi2 



(14) 



where z = e + iri is complex, and thus the imaginary part 
of the Green function is finite only inside the band, ±D, 
and is proportional to the semi-elliptic density of states. 
Thus the coupling constants are given by 



^TB 



(e) 



(15) 



In Sec. nil Bl we discuss the implementation of the po- 
larization, and explain that half-filled leads can be used, 
corresponding to e = 0. 



B. Modeling the polarization 

Full polarization of the leads is avoided for several rea- 
sons. Most prominently, full polarization decouples one 
spin species in a lead completely in the sense that the 
hopping matrix element between the lead and the reso- 
nant level is zero for all angles. Dealing with decoupled 
Hilbert spaces is undesired as it creates numerical prob- 
lems such as ill-conditioned matrices, making the numer- 
ical solution of the resolvent equations hardi^i^i^ 

Furthermore, there are single points where the model 
itself is ill-defined for full polarization. At the angles 
(f) — and (j) — TT the spin-flip process of the dot is in- 
active because of the prefactor sin (p. Due to the full 
polarization also the hopping matrix element for the mi- 
nority spin connecting the lead and the dot is zero. Thus 
the minority spin level is completely decoupled and hence 
has a constant occupation. The occupation of the major- 
ity spin level depends, however, on the occupation of the 
minority spin level through the on-site Coulomb interac- 
tion term Un-\ni. That is the properties of the model for 
these specific angles depend on the initial occupation of 
the minority spin level, and no unique stationary state 
exists. 

It should be noted that the qualitative behavior for 
(large) partial and full polarization are similar except 
for the problem for specific angles described above. It is, 
however, clear that a decreased polarization in the leads 
tends to wash out the spin dependence in the model, and 
in the limit of unpolarized leads all spin characteristics 
are lost. 

There is a certain freedom of choice in the modeling of 
the polarization. Although the polarization is a property 
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FIG. 5: (Color online) Sketch of the DMRG setup for the 
FAB model. Notice the implementation of the polarization 
through the hopping matrix elements, where ta,a and the on- 
site Coulomb interaction U are indicated in the figure. In the 
DMRG evaluation a single lead mapping is used in combi- 
nation with a momentum-space representation of the single 
tight-binding lead.^^ Here tk indicates a discretization depen- 
dent hopping to all states in the momentum-space. 



of the leads, it can be modeled by spin-dependent hop- 
ping matrix elements connecting the dot to the leadsj^ 
There are different approaches to modeling the polariza- 
tion of the leads and we have chosen the simplest one to 
implement in the DMRG setup. Rather than using spin- 
dependent filling in the leads, we use half-filled leads for 
both spin species, and model the polarization by modify- 
ing the hopping matrix elements connecting the leads and 
the dot. This is indicated in Fig. O where we show the 
DMRG setup using a momentum-space representation of 
the leads. This choice for the polarization simplifies the 
DMRG setup significantly as identical discretizations can 
be used for the two spin species in each lead such that 
the spin species are again treated equally apart from the 
polarization dependent hopping matrix elements. 



= t" 



(16) 



In all calculations presented, we use identical polariza- 
tions of the two leads, Pl — Pr = P, such that the 
coupling to the leads are identical when t'j^ = t^- 

In the remainder of the article we measure all ener- 
gies in units of the sum of the coupling constants at the 
equilibrium chemical potential, e = 0, 



(17) 



For the tight-binding chains this corresponds to mea- 
suring all energies in units of 
Eqs. (USD,®. 



C. Calculation of the conductance 

In order to obtain the conductance we make use of the 
Meir-Wingreen formula^'^ rather than the Kubo formula 
used in previous work.^^'^^ The evaluation of the spectral 
function at zero bias and for proportional couplings^^ can 
be done within a single lead framework, effectively halv- 
ing the lead size. For parallel magnetizations of the leads, 
the FAB model falls in the linear regime within this class 
of models, and thus the finite size scaling for the evalu- 
ation of the spectral function is significantly better than 
the evaluation of the Kubo formula, enabling faster and 
more accurate calculations)^ Using DMRG we thus eval- 
uate the two spin components of the full spectral function 
in separate calculations, and therefore need to recombine 
the spin resolved spectral functions into the total con- 
ductance. 



G(Co,0) = 



-T 



rA,(eo,0,c^ = O), (18) 



where the polarization enters through the hopping matrix 
elements ta,o-, and denotes the spin resolved spectral 
function of the dot. In this paper, we make the assump- 
tion that the hopping matrix elements between the leads 



and the dot are identical for both leads, t\ 



i^j, and that 
Pr. 



the polarizations in both leads are identical, Pl 

such that tL.a = tR^a = ta- 

In order to achieve the necessary precision in the 
DMRG calculations a momentum-space representation of 
the leads is used. Although the physics takes place at the 
Fermi level also energies well away from the Fermi level 
need to be represented properly, and the discretization 
scheme used should support this. We use a discretization 
of the momentum part of each lead consisting of a log- 
arithmic discretization that covers a large energy span, 
and switch to a linear discretization on the low-energy 
scale close to the Fermi leveL-'^ All the DMRG calcula- 
tions presented in this paper were performed using 55 
sites in the lead description, corresponding to 35 sites 
scaled logarithmically and 20 sites scaled linearly around 
the Fermi level 1^ 

It should be noted that by virtue of the DMRG method 
all interactions are rigorously taken into account. The 
approximation in the method presented lie in the use of 
a finite sized lead which can be benchmarked in the non- 
interacting limit, and as such the method used contains 
only controllable approximations. 



D. Results 

Using the momentum-space representation of the leads 
in the DMRG setup we have calculated the spectral 
function, and using the Meir-Wingreen formula in 
Eq. (fT5)) evaluated the conductance. For different values 
of the magnetic field strength B and the on-site Coulomb 
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interaction strength [/, we have calculated the conduc- 
tance versus the angle </> between the magnetic field and 
the polarization direction. In all examples, we keep the 
polarization in the two leads identical, Pl = Pb. — 0.8. 
As evident from the Hamiltonian, the model is sym- 
metric around (/) — n since cos(27r ~ cf)) — cos (j) and 
sin(27r — (j)) = — sin0 such that only the spin-flip term 
in H]j acquires an insignificant phase [see Eq. ([5])]. 
Therefore we confine our studies to angles in the interval 
(p S [0,7r], and the interval <p G [tt, 27r] is found by 
reflecting the results around (p = n. 

In order to determine the discretization needed for the 
leads, exact diagonalization calculations for the spectral 
function have been performed and, using Eq. ()18|) . com- 
pared to the NEGF results in the non-interacting limit, 
U ~ (not shown) i^i^ By virtue of the exact diagonal- 
ization, the only error present in this approach is the error 
due to the finite size of the leads. The results show ex- 
cellent agreement between the exact diagonalization and 
the Green function results for a range of parameter val- 
ues, and we conclude that the modeling of the leads is 
sufficient for resolving the model. 

Having benchmarked the DMRG setup in the known 
limit of [/ = we turn to the interesting regime of finite 
interactions. In Fig. [6] we show the results of the DMRG 
calculations ('-I-' symbols), keeping the bare level reso- 
nant, ^0 = 0. The calculations presented in each figure 
were performed keeping the strength of the magnetic field 
B fixed and varying the interaction strength U and the 
angle 0, where the specific parameter values are given in 
the plots. 

For the parameter regimes considered here the numer- 
ical zero-temperature DMRG results confirm the simple 
physical picture sketched in the introduction. That is 
the linear conductance versus the angle (j) shows anti- 
resonances for (/) = 7r/2 in the non-interacting limit, and 
a spin-valve behavior for strong on-site interaction. 

When S/r" decreases the maximum conductance in- 
creases as the levels move closer to resonance, but the 
qualitative behavior is the same when U is varied. Previ- 
ous work showed that for fully polarized leads and U — 0, 
the anti-resonances become sharper for decreasing B/T^ 
(see the left panel in Fig. 2 in Ref . Issl) . Due to the finite 
polarization, the sharpening of the anti-resonances is not 
very clear in the DMRG results. 

Finally we note that the cotunncling expression derived 
under the assumption i? ^ F", see Eq. (|Aip . reproduces 
the DMRG results fairly well already for B = 2r° (not 
shown) . 



IV. LINEAR RESPONSE: MEAN-FIELD 
SOLUTION 

Mean-field solutions are often problematic when con- 
sidering systems with only a few degrees of freedom such 
as, e.g., transport through quantum dots with only a few 
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FIG. 6: (Color online) The linear conductance vs. the an- 
gle for B = r°/2 and B = P". The symbols are 
the DMRG results and the full lines are the non-equilibrium 
Green function results using the Hartree-Fock approximation. 
We apply an elliptic density of states, using Eq. ((15} for the 
coupling constant with D = 2r°. The other parameters are 
P = Pl ~ Pr ~ 0.8 and = 0. The Green function result is 
exact in the non-interacting limit, U = 0, and is thus a rigor- 
ous benchmark for DMRG in this limit. The good agreement 
demonstrates the accuracy of the discretization, and confirms 
the capability of the DMRG. Furthermore there is also a sur- 
prisingly good agreement between these two methods for finite 
interactions as discussed in detail in Sec. IIVI 



levels contributing to the transport. In this case, the 
mean-field solution fails to describe Coulomb blockade^ 
and can lead to unphysical bistabilities at finite bias due 
to a sudden switch between different transport modes. 

Surprisingly, Fig. [S] shows that the mean-field version 
of the Hamiltonian from Eq. ([3]) actually reproduces the 
linear conductance results of the previous section rather 
well, which we will explain below. 

The mean-field version of the Hamiltonian in Eq. ([3]) is 
obtained by re-writing the interaction term as consisting 
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of a Hartree and a Fock-term 



Fock 



= U{dtd~^{dtd~^) + d\di{dtd^^)), 
= -U{d\di{d\d~^) + dtd^^idtdi)), 



(19) 



i.e., the replacement is Un^^n-^ -ffHartroc + ^Fock- Keep- 
ing only the Hartree-term gives a basis-dependent Hamil- 
tonian, and for the FAB-model the two spins are corre- 
lated (for 7^ 0, tt) due to the coupling to leads giving a 
non-vanishing Fock term.^^ 

From the mean-field Hamiltonian the linear conduc- 
tance can be obtained using the non-equilibrium Green 
function formalism, and the calculation is similar to the 
non-interacting calculation in Ref. Issl . However here the 
non-interacting dot Green function is replaced with the 
Hartree-Fock dot Green function 



u;-iCo-B + U{dldi)) 
+U{dUi) 



+u{dl 

(Co + S ^ 



-U{d\d~^)) 



(20) 



The generalized occupations have to be calculated self- 
consistently through the relation 



with the lesser Green function being 



(21) 



G<(c^) = i,fL{Lo)G^{uj)TLG''{uj)+xfB{uj)G'^{Lo)TRG''{Lo), 

(22) 

where the expressions for the coupling constants, Fq, are 
derived in App. [B) 

The zero-temperature linear conductance is obtained 



^47.55 



G 



MP 



h 



Tt[G\Q)TlG'-{Q)Tr]. 



(23) 



The results for B — r°/2 and i? = F" are shown in 
Fig. [6] together with the DMRG results. It is clearly seen 
that for S = r° the results agree almost exactly, but for 
B < r° deviations start to appear, especially for angles 
around 7r/2. This rather surprising success of the mean- 
field solutions can be understood by an inspection of the 
occupations, see Fig. [T] 

First we notice that the couplings to the levels de- 
pend on the angle see Eq. ([10]), and consider the non- 
interacting case U = Q. At = the spin-t level is 
roughly 80% occupied due to a large broadening, but 
the spin-| level is almost unoccupied due to the large 
^/^L II + II) ~ At (/) = TT the spin-f-level is al- 
most decoupled and the occupation goes to 1, whereas 
the occupation of the spin-| level is ~ 20%. For the in- 
termediate angles, there is a smooth cross-over between 
the two regimes, see the black curves in Fig. [T) When 
calculating the occupations self-consistently for finite U , 




7r/2 
Angle 



FIG. 7; (Color online) The mean-field occupations vs. (j> for 
B — r''/2 and four different values of the on-site Coulomb 
interaction. The upper bunch of curves are for (n;^) and the 
lower ones are for (nj). The rest of the parameters are as in 
Fig. El 



the overall trend does not change: (n^) remains almost 
identical, and (n^) decreases for increasing [/, see Fig.[71 
So the reason why the mean-field solutions performs 
relatively well (at least at zero temperature) is that 
when one of the levels fluctuates the most, the other 
level has an occupation being either approximately zero 
or 1, i.e., the products of the fluctuations vanishes. 
For this specific model, the success of the mean-field 
solution is due to the combination of the split levels and 
the angle dependent couplings, which means that the 
term (n^ — (n|))(nj- — {n^)), neglected in the mean-field 
Hamiltonian, remains small for all angles. The largest 
deviation between the mean-field results and the DMRG 
results is expected for the angles around = 7r/2, which 
is also observed in Fig. [SI 

At elevated temperatures or for even smaller mag- 
netic fields, the mean- field solution is expected to per- 
form worse due to larger fluctuations of the occupa- 
tions, but for finite temperatures no exact results are cur- 
rently available for comparison. Furthermore, while we 
have shown here that the mean-field solutions are some- 
what fortuitously reliable to describe linear response, we 
demonstrate in Sec. lVIII that they fail for the FAB model 
for finite bias. 



V. LINEAR RESPONSE: HUBBARD-I 
APPROXIMATION 

Another popular and widely used approximation is the 
Hubbard-I approximation (HIA) which corresponds to a 
decoupling of equations of motion for the Green's func- 
tions at a higher level of the hierarchy than the simple 
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mean-field approximation used in the previous section. 
Therefore, one could expect better performance of this 
approximation compared to the mean-field approxima- 
tion. This expectation is further supported by the fact 
that HIA naturally arises as the lowest-order approxima- 
tion in the theories developing perturbation series around 
the state of the isolated system, i.e. around the so-called 
"atomic limit" . Despite the lack of existence of the stan- 
dard Wick theorem due to the non-Gaussian nature of the 
unperturbed system with arbitrary correlations, system- 
atic (even renormalized) perturbation theory reportedly 
exists^^i^ and HIA is in a certain sense the lowest order 
in that expansion. Indeed, it turns out that HIA can 
rather simply and yet correctly describe the non-trivial 
effects of single-level broadening in the Coulomb block- 
ade regime of transport studied by other sophisticated 
methods^i^ and confirmed experimentally.— Therefore 
HIA appears to be the solution to the problem of a simul- 
taneous description of interference/broadening and inter- 
action. 

However, this optimistic picture breaks down as soon 
as more complicated systems are addressed, namely, 
any system where coherence between different transport 
channels is involved. The FAB model is one such exam- 
ple and we will demonstrate the breakdown of HIA for 
this model. A similar situation arises for the problem of 
a single spin-degenerate level with local Coulomb inter- 
action coupled to superconducting leads in the so-called 
TT-junction regime, where the Joscphson supercurrent is 
observed. The HIA is known to fail to predict this the- 
oretically and experimentally well-established fact, see 
e.g. Ref. [31I, Fig. 4a. The situation in the case of the 
FAB model is even more severe since the HIA failure is 
not just physical, as in the superconducting case, but the 
results are even mathematically inconsistent. Fundamen- 
tal analytical identities such as (C)^ = G° and, conse- 
quently, the hermiticity of the density matrix is broken 
within the HIA. 

To demonstrate this explicitly, we use the Hamilto- 
nian in Eq. Q and perform analogous derivations as in 
Refs. 26s ^ for simpler models without any magnetiza- 
tion at all or without the local splitting, respectively. We 
arrive at the matrix equation 



G(£)=g(e) + g(e)S(£)G(£) 



(24) 



for the causal Green function G^a' (e) (in energy repre- 
sentation) of the central dot. Here 



1 



(25) 



ka 



is the self-energy matrix due to the coupling to the con- 
tacts, and we introduced the auxiliary Green function 



I £-^^-C/(l-(rftdj)) 

-U{dU~^) e-e.j-U{l-{dU-^)) I ' ^^^1 



where ^^^j- = =p i? are the level energies in the basis 
parallel to B as given in Eq. ([6]). For vanishing split- 
ting, B — Q, these are precisely the Eqs. (15)-(21) of 
Ref. (2^. The retarded and advanced Green function are 
then again obtained by G'^l'^ie) = G(e ± iO+), respec- 
tively. In addition, we define the retarded and advanced 
components of the other functions in the same way. 

A general condition relating the retarded and advanced 
Green's functions is (G"")^ = G". As Eq. ([M]) implies 



Qr/a 



(27) 



and (S'')'i' = T." from Eq. ^ this requires (g'')''' = g". 
This is, however, only compatible with Eq. ([26|) if the 
off-diagonal elements are mutually complex conjugated. 
For non-zero coherences (fi~c?|) = (djiij^)* 7^ and finite 
level splitting ^rj: ^ pertinent to the FAB model this 
condition is not satisfied and, thus, HIA breaks the nec- 
essary mathematical condition for the Green's functions. 
Similarly, self-consistent evaluation of the (generalized) 
populations in the spirit of Eq. (j2ip in thermal equilib- 
rium (since we study currently the linear response only) 
would yield a non-hermitian density matrix, yet another 
mathematical problem stemming from the inconsistency 
of the HIA equations. 

It should be stressed that the inconsistency only shows 
up if both the local level splitting and the coherences be- 
tween different spin states (stemming from non-collinear 
magnetization arrangement) are present. Therefore, the 
inconsistency has apparently not been noticed before^i^ 
since previously employed models do not contain both 
necessary ingredients. 

Nevertheless, the FAB model is both mathematically 
and physically a realistic model and the failure of HIA 
reveals problems inherent in that approximation. As 
mentioned before, the HIA is in some sense the lowest- 
order expansion in reportedly systematic theories based 
on Hubbard operator Green's functions. '^^ One can show 
that the same problems carry over to the many-body for- 
malism, see Ref. for a pedagogical overview of the con- 
nection between the formulation of HIA in the standard 
as well as the many-body formalism. 



VI. FINITE BIAS: 2VN APPROACH 

So far we have only considered linear response and 
zero temperature. In this section we apply a density ma- 
trix formalism developed in Ref. ;46; , giving access to the 
regime of finite bias and finite temperature. The method 
works in a basis of many-particle eigenstates for the dot 
Hamiltonian, thereby including all interactions on the 
dot exactly. Correlated transitions between the lead and 
the dot states with up to two different electron states 
are included exactly, which suggest the notation second 
order von Neumann approach (2vN). By solving the re- 
sulting set of equations for the steady-state, a certain 
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class of higher-order processes are also included. Inter- 
ference effects are also included by a full treatment of the 
nondiagonal density matrix elements. 

The dot part of the Hamiltonian, Hd, has four 
many-particle eigenstates {|0), |t): |2)}, where |2) = 
|0) , and the energies are 0, = ^q — B,E~^ = + -B 
and E2 = E~ + E~ + U, respectively. Inserting a com- 
plete set of dot states in the tunneling Hamiltonian from 
Eq. ([8]) gives 

Ht= J2 {na{kaa)\b){a\co..k,+cl,Ja){b\T:,{kaa)) 

k(Ta.,ab 

(28) 

where a, b denotes the dot many-particle eigenstates and 

Tba{kaa) = Y.^.^aM^'^t'^^^^U"''^ couplings be- 

tween these states and the lead states. 

Inserting these coupling matrix elements and eigenen- 
ergies in Eqs. (10), (11) in Ref. [H gives a closed set of 
equations for the elements of the reduced density ma- 
trix, which can be solved numerically. As for the DMRG 
calculation we assume the leads to be 80% polarized. Sig- 
nificantly larger polarizations (« 1) are difficult to handle 
numerically, especially in the linear conductance regime. 
This also limits the values for the bias voltages and tem- 
peratures which can be used. 

For the finite bias calculations presented here, we use 
a constant density of states, so that effects due to the 
change of the chemical potentials are not superimposed 
by changes in the contact couplings. For numerical pur- 
poses, we implement a finite band width with elliptically 
shaped edges at 0.951? < |e| < D, where D surpasses all 
other relevant energy scales. In all plots, the bias voltage 
V is applied symmetrically, — — /^fl — eV/2. 

Before considering finite bias, we have carefully in- 
spected the results for low bias {eV/T'^ — 0.05 — 0.1) 
and low temperature {ksT/T^ = 0.05 - 0.1) for B/r° = 
0.5,1,2 for all angles (not shown). For the non- 
interacting case, the exact NEGF results are reproduced 
for all parameters tested. From the low-bias results we 
have extracted a numerical value for the linear conduc- 
tance, and the results show almost quantitative agree- 
ment with the exact DMRG results for all tested values 
of U and all angles. The discrepancies can be attributed 
to the (small but) finite bias and the finite temperature 
used in the density matrix calculation. We conclude that 
the 2vN approach is capable of describing the effect of 
the interactions and the coherence in the low-bias regime 
for the model system considered. 

Fig. [8^) shows the current versus bias voltage for differ- 
ent angles, where the bare level is on resonance, = 0, 
B ~ 2r° and an on-site Coulomb interaction U — 8r". 
Shoulders in the current are expected if half the bias 
matches the single-electron transition energies in the dot. 
This happens at eF/F" = 4 for the transitions |0) |t) 
and |0) ^ ll), at eV/V" = 12 for |I) ^ |2), and finally 
at eV/r" = 20 for |f) -> |2). 



In the low-bias regime, eV/T^ < 4, the current is sup- 
pressed when the angle (j) is increased from to tt due to 
the spin- valve effect (see also, e.g., Fig.|6]for large U /F°). 

In the intermediate regime 12 < eV/T'^ < 20, the cur- 
rent shows a very pronounced angular dependence, with 
a significant current drop between ip — and (p — tt. For 
(/) « TT one even detects negative differential conductance 
around eV — 12r°. In this region, the lower spin level "f 
can be filled not only by the process |0) — > |t) but also 
by li) ^ |2). It is thus more likely to be filled com- 
pared to lower biases eV < 12r°. In the case </> = 0, t is 
aligned with the lead polarization and therefore exhibits 
high tunneling rates. Thus its increased occupation prob- 
ability goes along with an increase of current around 
eV = 12r°. In contrast, for = tt, "f is pointing against 
the lead polarization and therefore has a low tunneling 
rate, explaining the drop of current around eV = 12r*'. 
For intermediate angles there is a smooth cross-over be- 
tween the two limits. Here the non-diagonal elements of 
the density matrix are non-vanishing and quantum co- 
herence plays a role, as the two dot states are superpo- 
sitions of the lead spins. Off-diagonal elements are also 
important to include in transport through dots coupled 
to non-collinear ferromagnetic leads, even in the absence 
of an applied magnetic field^^i^ 



VII. FINITE BIAS: MEAN-FIELD NEGF AND 
MASTER EQUATION 

Figure [8b) shows the current versus applied bias cal- 
culated using the Hartree-Fock mean-field version of the 
Hamiltonian, see Eqs. (|19p - ([22|) and App. [Bl within the 
NEGF formalism using a self-consistent calculation of the 
occupations, see, e.g., Ref. WA The parameters are iden- 
tical to Fig. [8^) in order to allow for a direct comparison. 

In the low-bias regime {eV < 4r"), where the aver- 
age occupations (n^), (nj) are close to one and zero, re- 
spectively, the results agree with the 2vN formalism of 
Sec. I VII This goes well with the observation from sec- 
tion |TVl that the conductance is well reproduced within 
the mean-field model. 

In the region 4 < eV/F" < 12, we find {n:^) « 0.5 as 

the lower energy state t is in the window between the 
Fermi levels with symmetric coupling to both contacts. 
Thus the higher energy state | has the energy B + U/2 = 
6r° which is above the emitter Fermi level and {uj) « 

0. If f is aligned with the lead polarization (i.e. (p — 
0), the current is larger, while it is low for (f> — t: as 
seen in Fig. [SJd). Finally, at eV/T^ > 12, the upper 
level i enters the window between the Fermi levels, takes 
an average occupation (nj) ~ 0.5, and contributes also 
to the current. The repulsion of the levels by U/2 is 
an artefact of the mean-field model, and correspondingly 
neither the current values nor the shoulders agree with 
the more detailed 2vN results shown in Fig. [5^). 

Fig. [5J;) shows the corresponding result from the mas- 
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FIG. 8; (Color online) The current versus bias voltage for 
five different angles <j) obtained using the 2vN density matrix 
formalism (a), the mean-field Hartree-Fock approach (b), and 
the master equation approach (c). The parameters are: = 
0, B = 2r°, U = 8T°, ksT = 0.ir°, Pl= Pr = 0.8 and the 
bias is applied symmetrically. For the 2vN-method we used a 
constant density of states with a half band- width D — 20r'', 
see the main text, while for the Hartree-Fock calculation the 
wide-band limit is applied, i.e. the real-part of the self-energy 
was neglected. 



ter equation approac h^^'^^i^^ where the occupations of 
the many-particle states are determined by electron hop- 
ping processes to and from the leads. While this ap- 
proach does not provide any current for low biases, where 
cotunncling dominates the transport, it provides reliable 
results for the current plateaus. In particular, the oc- 
currence of negative differential conductance for « tt 
for biases around 12r°/e is confirmed. Note that the 
presence of pronounced steps of width fc^T is due to the 
entire neglect of broadening in this approach. Similar 
results are obtained by taking into account nondiagonal 
density matrices within the IvN approach^ (not shown) . 

VIII. CONCLUSION 

In this article, we provided a full description of the Fer- 
romagnetic Anderson model with applied magnetic field 
B (FAB). We have successfully implemented the density 
matrix renormalization group (DMRG) method, which 
provides the linear conductance for arbitrary strength of 
the on-site Coulomb interaction and arbitrary level split- 
ting. The data interpolate between the known results of 
non-equilibrium Green functions (NEGF) for zero inter- 
action and the cotunncling results for large level splitting. 
A key result is the strong suppression of conductance 
with increasing on-site Coulomb interaction if the mag- 
netic field on the dot is opposite to the lead polarization. 

The DMRG results can serve as a benchmark for dif- 
ferent approaches, where we find that both the second 
order von Neumann (2vN) approach and the NEGF ap- 
proach with mean-field interaction give reliable results for 
the conductance. While the 2vN approach also provides 
plausible results for finite bias, the mean- field NEGF fails 
due to the wrong treatment of partially occupied states. 

For finite bias the 2vN approach predicts a strong de- 
pendence of the current on the direction of the magnetic 
field in the intermediate bias region. Here negative dif- 
ferential conductance is predicted if the magnetic field 
on the dot is opposite to the lead polarization. This fea- 
ture can also be qualitatively obtained from the simpler 
master equation approach. 

Finally we have shown that the Hubbard I approxima- 
tion leads to unphysical results for this particular model. 
This shows that the FAB model constitutes a sensitive 
test case for different approaches due to its involved inter- 
play between interference, broadening, and interaction. 
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APPENDIX A: COTUNNELING EXPRESSION 

The expression for the cotunneling current in Eq. ([2|) 
can easily be generalized to arbitrary polarization and 
off- resonant transport, $o 7^ 0. For identical polarizations 
of the leads, Pl — Pr = P, it reads 



G 



„2p0 pO 
cotun L'^ L 



2T:h 



l + P 



cos2(0/2) sin2(0/2) 



+ 2 



(1 + P)(1-P) 



Co - S ^o + B + U 



- sin((/«/2) cos{(t>/2) sin(0/2) cos((/«/2) 



^o + B + U 
I 



(0/2) , cos2(0/2) 



-B 



^0 



B + U 

(Al) 



APPENDIX B: EQUATIONS FOR THE NEGF 
SOLUTION 

For completeness we present here the equations for 
the non-equilibrium Green function calculations within 
the Hartree-Fock approximation. Using the equation-of- 
motion technique, the Green functions in the diagonal 
basis are^ 

^■-(e) = (Gh^ - S2'^(e) - ^Ti^)y' , (Bla) 
- T.i^^)^'^\i^M\'9ZkAe)n,a' , (Bib) 

kCT 

where = (e - Lm ^ *0+)-i and r„ = 

—i [S^ — SJ^]. The Green function G^p(e) is stated in 
Eq. 

For the tight- binding chain of Sec. IIIII with elliptic 
bands we find 



r„(e) =2^^|ia,fe.|'<5(e-ea,fea) 

kcr 



(B2) 



where Eqs. have been used. The ff -component 

of the self-energy becomes 



K"(s)]ff =E 



1, f de' 

--il + P^cos^)J - 



e-e'± iO+ 



(B3) 



where we used the definition of the polarization [see be- 
low Eq. (O]. 

The principal part of the integral can be found analyt- 
ically for the elliptic density of states 



-D 



2tt e-e' 



2D 



for lei < D, 



(B4) 



where it has been used that rQ(0) = 4|i°p/I?. 

The other components are calculated similarly, and the 
full expression for the self-energy becomes 



1 / 1 + Pa COS ( 

Pa sin (j) 



Pa sin (j) 

I ~ Pa COS < 



£r„(0) iTaie) 
2D ^ 2 



(B5) 



After a self-consistent evaluation of the Green func- 
tions, the current can within the mean-field approxima- 
tion be evaluated as 



1 /""^ 

2TThJ_^ 

X [/L(e) -/«(£)], 



(B6) 



with fa{e) = 1/ [e(^-t^^)/kBT + i]^ Q, 
Fermi function. 



L, i?, being the 
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